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DEC's System Research Center is a new group, still recruiting its initial members and laying plans for 
long-term work. The business and technology objectives of DEC require a strong commitment to 
research. We join two other corporate research groups in meeting that commitment. 

SRC's role is to design, build, and use new digital systems five to ten years before they become 
commonplace. Our purpose is to advance both the state of knowledge and the state of the art. 

SRC will create and use real systems in order to investigate their properties. Interesting systems are too 
complex to be evaluated purely in the abstract. Our strategy is to build prototypes, use them as daily 
tools, and feed the experience back into the design of better tools and more relevant theories. Most of 
the major advances in information systems have come through this strategy, including time-sharing, the 
Arpanet, and distributed personal computing. 

Among the areas in which SRC will be building prototypes during the next several years are applications 
of high-performance personal computing, distributed computing, communications, databases, 
programming environments, system-building tools, design automation, specification technology, and 
tightly coupled multiprocessors. 

We will also do work of a more formal and mathematical flavor; some of us will be constructing theories, 
developing algorithms, and proving theorems as well as designing systems and writing programs. Some 
of SRC's theory work will be in established fields of theoretical computer science, such as the analysis of 
algorithms, computational geometry, and logics of programming. In other cases, we expect to explore 
new ground motivated by problems that arise in our systems research. 

DEC has a commitment to open research. The improved understanding that comes with widespread 
exposure seems more valuable than any transient competitive advantage. SRC will freely report results 
at conferences and in professional journals. We will encourage visits by university researchers and 
conduct collaborative research. We will actively seek users for our prototype systems. To facilitate 
interchange, we will develop systems that run on hardware available to universities and work out ways of 
making our software available for academic use. 
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Capsule review: 



When the user clicks a mouse button on a multi-window 
graphics workstation, the system must determine which window 
the cursor was in. This is one instance of the point 
location problem: Given a partition of the plane into regions 
and a sample point, identify the particular region in which 
the point lies. The point location problem has received lots 
of attention from theoreticians, and many algorithms have 
been proposed. The new algorithm presented in this paper is 
the third to win the theory prize, that is, to achieve the 
best possible bounds on time and space performance. Of those 
three, it is the first that a programmer might reasonably 
adopt for use in a real system. 

The new algorithm is best thought of as a refinement of the 
previous practical champ, due to Lee and Preparata. That 
algorithm threw away information at various points. The 
authors of this paper discovered an ingenious way to preserve 
essentially all of this information while increasing the 
storage requirements only moderately. 

The authors pay attention to the practical details of data 
structure design and to the size of constant factors. But 
they keep the discussion at a fairly high level — don't 
expect complete Pascal programs. The geometric ideas 
involved are clearly portrayed by the many pretty figures. 
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Abstract: Point location, often known in graphics as "hit detection," is one of 
the fundamental problems of computational geometry. In a point location query we 
want to identify which of a given collection of geometric objects contains a particular 
point. Let S denote a subdivision of the Euclidean plane into monotone regions by 
a straight-line graph of m edges. In this paper we exhibit a substantial refinement 
of the technique of Lee and Preparata for locating a point in S based on separating 
chains. The new data structure, called a layered dag, can be built in O(m) time, uses 
O(m) storage, and makes possible point location in O(logm) time. Unlike previous 
structures that attain these optimal bounds, the layered dag can be implemented 
in a simple and practical way, and is extensible to subdivisions with edges more 
general than straight-line segments. 

Keywords: point location, monotone polygons, planar graphs, partial order, graph 
traversal, computational geometry. 

1. Introduction 

Point location is one of the fundamental problems in computational geometry. In the two- 
dimensional case, we are given a subdivision S of the plane into two or more regions, and 
then asked to determine which of those regions contains a given query point p. If the same 
subdivision is to be used for a large number of queries, as is often the case, we can reduce 
the total cost of this task by preprocessing the subdivision into a data structure suitable for 
the search. We will measure the performance of a proposed solution to this problem by three 
quantities, the preprocessing cost P, the storage cost S, and the query cost Q. 

Optimal solutions for this search problem have been known since Lipton and Tarjan [LT] and 
Kirkpatrick [Ki]. For a subdivision S with m edges these optimal solutions simultaneously 
attain S = O(m), Q = O(logm) and, under certain assumptions, also P = 0(m). The 
Lipton- Tarjan method is based on their graph separator theorem, and so far has been only of 
theoretical interest. Kirkpatrick 's method, which consists of building a hierarchy of coarser 
and coarser subdivisions, is implementable, but still the implied constants are so large as 
to make current implementations of little practical interest. In addition, neither of these 
techniques seems to extend in a natural way to curved-edge subdivisions. 

The work of Herbert Edelsbrunner was supported by the Austrian Fonds sur Forderung der wissenschaftlichen 
Forschung. Jorge Stolfi, who is on leave from the University of Sao Paulo, was partially supported by a grant 
from Conselho Nacional de Desenvolvimento Cientffico e Tecnol6gico (CNPq). Part of this work was done 
while the last two authors were at the XEROX Palo Alto Research Center. 
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Historically, Dobkin and Lipton [DL] were the first to achieve O(logm) query time, while us- 
ing 0(m 2 ) space. Their method was subsequently refined by Preparata [P] so that 0(m log m) 
space suffices. Later Bilardi and Preparata [BP] again gave a refinement that achieves O(m) 
space in the expected case, while retaining 0(m log m) space and O(logm) query time in the 
worst case. These solutions are applicable to curved-edge subdivisions and seem to admit of 
efficient implementations. 

A substantially different approach was taken by Shamos [S] and led to the well-known point 
location paper of Lee and Preparata [LP], based on the construction of separating chains. 
This data structure attains P = 0(m log m), 5 = O(m), and Q — 0(log 2 m). The constants 
of proportionality in these expressions are quite small, and the algorithms, particularly the 
query one, are simpler than those of Kirkpatrick. While Kirkpatrick's algorithm requires the 
regions to be triangular, that of Lee and Preparata works for regions of a more general shape 
(monotone polygons, which include the convex ones). Recently Chazelle [C2] described a 
variant of this solution as a special instance of a general method for "searching in history" . 
His structure needs the same amount of space and time. These techniques again are applicable 
to curved-edge subdivisions, although this possibility was not examined. 

In a separate development, Edelsbrunner and Maurer [EM] came up with a space-optimal 
solution that works for general subdivisions, and even for sets of arbitrary non-overlapping 
regions. The query time is Q = 0(log 3 m), but it can be improved to Q = 0(log 2 m) for 
rectangular subdivisions, where the structure becomes especially simple. For this reason a 
generalization to rectangular point location in higher dimensions has also succeeded; see 
Edelsbrunner, Haring, and Hilbert [EH]. 

The purpose of this paper is to show an elegant modification to the separating chain method 
of Lee and Preparata that yields a new optimal point location algorithm for monotone sub- 
divisions of the plane. The algorithm is based on a new data structure called the layered dag. 
In this new structure the separating chains built into a binary tree by Lee and Preparata are 
refined so that (1) once a point has been discriminated against a chain, it can be discrimi- 
nated against a child of that chain with constant extra effort, and (2) the overall storage only 
doubles. The layered dag simultaneously attains S — 0{m) and Q = 0(log m). An additional 
insight allows us to build this dag (as well as the original Lee-Preparata tree) in only O(m) 
time. Not only is this new structure optimal with respect to all of preprocessing, space, and 
query costs, but in fact a simple implementation, with small constants of proportionality, is 
possible. Like its Lee-Preparata predecessor, it also extends to curved-edge subdivisions. 

In the organization of the paper we have adopted the policy that each new data structure is 
first introduced by how it is to be used, and then by how it is to be constructed. We have 
placed emphasis throughout on implementation considerations, as well as the underlying the- 
ory. Section 2 describes the basic notions surrounding monotone polygons and subdivisions. 
Section 3 shows how non-monotone subdivisions can be made monotone. Sections 4, 5, and 
6 introduce a partial ordering of the regions and its use in getting a complete family of 
separators. Section 7 reviews the Lee-Preparata structure, while sections 8 and 9 introduce 
the layered dag, and explain its use in point location and its construction. Section 10 gives 
an algorithm for constructing a complete family of separators based on a traversal of the 
subdivision. Two implementation issues are taken up in sections 11 and 12; these may be 
omitted on a first reading. Section 11 describes some bit-twiddling trickery used to give us 
linear preprocessing time, while section 12 disease how subdivisions can be traversed without <*\scuss*-3 
auxiliary storage. Finally section 13 contains some further applications of the layered dag to 
problems in computational geometry. 
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2. Monotone polygons and subdivisions 

An interval is a convex subset of a straight line, i.e., the whole line, a ray, a segment, a single 
point, or the empty set. An interval is proper if it contains more than one point, and is open 
if it does not contain its endpoints (if any). A subset of the plane is said to be monotone if 
its intersection with any line parallel to the y axis is a single interval (possibly empty). 

We define a polygon as an open, connected, and simply connected subset of the plane whose 
boundary can be partitioned into finitely many points (vertices) and open intervals (edges). 
Note that, according to these definitions, polygons may have infinite extent. A subdivision 
is a partition of the plane into a finite number of disjoint polygons (regions), edges, and 
vertices; these parts are collectively called the elements of the subdivision. It follows from 
the definitions that every edge is on the boundary between two regions (not necessarily 
distinct), that every vertex is incident to some edge, that every endpoint of an edge is a 
vertex, and that every region (unless there is only one) has some edge on its boundary. From 
these facts and from Euler's theorem, we can conclude that a subdivision with m edges has 
0(m) vertices and regions, thus justifying our use of m as the measure of problem size. 

A subdivision is said to be monotone if all its regions are monotone and it has no vertical 
edges. The last condition is a technical one: it is imposed only to simplify the proofs and 
algorithms and can be removed with some care. Figure 1 illustrates a monotone subdivision 
(arrowheads denote edges extending to infinity). 




Figure 1. A monotone subdivision. 

With minor caveats and modifications, monotone subdivisions include many interesting sub- 
cases, such as triangulations, subdivisions of the plane into convex pieces, and so forth. The 
lemma below shows that monotone subdivisions are precisely the "regular planar straight-line 
graphs" as defined by Lee and Preparata [LP]: 
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Lemma 1. A subdivision with no vertical edges is monotone if and only if every vertex is 
incident to at least two distinct edges, one at its left and one at its right. 

Proof: The "only if" part is easy, since if, for example, all edges incident to a vertex v 
pointed to the right, then the region to the left of v would have a disconnected intersection 
with a vertical line through v. 

For the converse, assume R is a region that is not monotone, and let / be a vertical line 
whose intersection with R consists of two or more components, as in figure 2. Since R is 
connected, any two components I\, I2 of the intersection are connected by a simple path 
jt in R. We can always choose I\, I2, and jt so that I\ is above I2, and jr does not meet 
/ except at its endpoints pi € A and p 2 G /j. Then it and the interval [pi,P2] delimit a 
closed and bounded region S of the plane. 




Figure 2. 

Suppose jt lies to the left of /. The boundary of R must enter S at the lower endpoint of 
I\ (and also at the upper endpoint of I2). Since the boundary cannot meet the path ir, 
there must be some vertex of the subdivision in the interior of S. Let v be the leftmost 
such vertex; clearly all edges incident to v lie to the right of it, and we have proved the 
lemma. A similar argument holds if jt lies to the right of /. I 

Therefore, we can check whether any subdivision is monotone by enumerating all the edges 
incident to each vertex, and checking whether they satisfy lemma 1. Each edge will be ex- 
amined at most twice, and each vertex once, so this algorithm runs in 0(m) time. 
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3. Regularization 

Lee and Preparata [LP] have shown that an arbitrary subdivision with m edges can be refined 
into a monotone subdivision having < 2m edges in 0(m log m) time, a process they termed 
regularization. (They define a vertex as being regular if it satisfies the condition of lemma 1.) 
For the sake of completeness we reproduce their algorithm here, in a slightly different setting. 

The task of the regularization process is to add new edges to the subdivision so as to make 
each vertex regular. Each new edge must connect two existing vertices (or extend from an 
existing vertex to infinity), and may not intersect any other edges. In other words, we can 
connect two vertices only if they are visible from each other. See figure 3 for an example. 



The regularization algorithm is based on the sweeping line paradigm that has been used 
extensively in computational geometry. We imagine that a vertical line sweeps the plane 
from left to right. We call active those vertices, edges, and regions currently intersecting the 
sweeping line. The list of active elements changes only when the line passes through a vertex, 
at which time some edges may end and others may start. Except at those moments, the 
active edges have an obvious vertical ordering, and cut the sweeping line into one or more 
active intervals. To each active interval we assign a generator, which is the last vertex swept 
over that happened to lie on or between the two active edges bounding that interval. This 
vertex may be a left endpoint of those edges, as figure 4 shows, or it may be an irregular 
vertex with no right-going edge, as figure 5 shows. 




Figure 3. Regularizing a non-monotone subdivision. 




generator 



I sweeping 
1 line 



Figure 4. 



Figure 5. 
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Lemma 2. The generator of an active interval is visible from any point in that interval. 

Proof: By definition, the hatched region in figures 4 and 5 is free of vertices and edges. I 

The regularization algorithm simulates this line sweep. We start by sorting all vertices of 
the subdivision by their x-coordinates. We will assume these coordinates are all distinct; if 
necessary, we can enforce this condition by rotating the subdivision through a sufficiently 
small angle. We also augment the subdivision with two dummy vertices at x = ±oo, to 
which are incident all edges having infinite left and right extent, respectively. During the 
sweep, we maintain a balanced search tree whose leaves represent the active intervals of the 
current sweep line (and the associated generators), ordered by y-coordinate. Conceptually, 
we imagine that the sweep line jumps from one corridor between successive vertices to the 
next, because within a corridor neither the ordering nor the generators of the active intervals 
can change. 

Consider what happens when we pass a vertex v with / left-going and r right-going edges. At 
that moment we must delete the / + 1 intervals bounded by edges that end at v, and insert 
in their place the r + 1 intervals bounded by edges that start at v. In particular, if / = 0, 
we must delete the currently active interval in which v lies, which is bounded by the edges 
immediately above and below v. Similarly, if r = 0, we must add one new interval bounded 
by those two edges. Updating the generators is straightforward: the vertex v simply becomes 
the generator of the r + 1 new intervals. 

When an active interval is about to be deleted, we check whether a new regularizing edge 
should be added to the subdivision. If the generator u of that interval has no right-going 
edges, or the new vertex v has no left-going ones, we add a new edge between u and v. See 
figure 6. By lemma 2, this edge does not intersect any others. 




Figure 6. Adding new edges. 

When the sweep line reaches the vertex at x = +co, all irregular vertices will have been fixed 
in this way. The number of new edges is at most the original number of irregular vertices, 
which in turn is at most m. The running time is dominated by the cost of sorting of the 
vertices and manipulating the balanced tree, which is O(mlogm) in the worst case. The 
problem of whether regularization can be done in time faster than O(mlogm) remains open. 
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4. The vertical ordering 

We denote by UA the orthogonal projection of a set A on the x-axis. The projection UR 
of a monotone polygon is an open interval of the x-axis, whose endpoints, if any, are the 
projections of at most two extremal vertices of R. The frontier of R is the boundary of R 
minus its extremal vertices. This frontier consists of two disjoint pieces, the top and the 
bottom of the region R. Each is either empty, or a connected polygonal line, which may 
extend to infinity in one or both directions. 

Given two subsets A and B of the plane, we say that A is above B (and write A » B) if 
for every pair of vertically aligned points (x, y a ) of A and (x, yb) of B we have y„ > y&, with 
strict inequality holding at least once. In this case we also say that B is below A (and write 
B « A). Some of these concepts are illustrated in figure 7. For the rest of this section, we 
will restrict » to the elements of a fixed monotone subdivision S, with n regions and m 
edges. Then » has several interesting properties listed in the lemmas below (straightforward 
proofs will be omitted). 




ITB 



Figure 7. A » 5, 

Lemma 3. For any two elements A and B of a monotone subdivision, A » B if and only 
if A ^ B and there is some vertical line segment with lower endpoint in B and upper 
endpoint in A. I 

Lemma 4. Any two elements A and B of a monotone subdivision satisfy exactly one of 
A = B, A « B, A » B, or IIA n I1B = <f>. I 

If three elements A, B and C are intercepted by a common vertical line, then from A « B 
and B « C we can conclude A « C. Therefore, the relation « is transitive (and in 
fact a linear order) when restricted to all the elements intercepted by a given vertical line. 
Transitivity may not hold without this restriction, but the following property will be true in 
any case: 
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Lemma 5. The relation « is acyclic. 

Proof: Suppose not. Let E 0 « E 1 « E 2 « • • • « E k — E 0 be a cycle of « with 
minimum length, where each Ei is an element of the subdivision S. From lemma 4, it 
follows immediately that k > 2. 

The x-projections of any two consecutive elements Ei and Ei +1 in this cycle must have 
some abscissa x ; in common. If for some t we had III^ c IT^-i, then the vertical line 
through Xi would meet Ei- X , E it and E i+1 ; transitivity would hold, and we would have 
« Ei+i. See figure 8. But then we could omit Ei and still get a cycle, contradicting 
the assumption that k is minimum. For analogous reasons, we cannot have TLEi C 1125,-+ 1. 




Figure 8. 

Let Ei be one of the "leftmost" cycle elements, i.e., such that none of the intervals UEj 
contains a point to the left of UEi. The projections n£,_i and IlEi+i of its neighbors must 
both meet UEi, and, by what we have just said, extend beyond its right endpoint. But then 
there would be again a vertical line meeting all three elements, implying Ei-i « Ei +i 
and k is not minimum. We conclude that no such cycle exists. I 

We say that a region A. is immediately above a region B (and write A >- B) if A » B and 
the frontiers of the two regions have at least one common edge; see figure 9. 




Figure 9. A >- B. 

In general, the >- relation is stronger than », but the following is easily seen to be true: 

Lemma 6. The transitive closure of » (restricted to the regions of a subdivision) is the 
same as that of >-. I 
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5. Separators 

A separator for a subdivision S is a polygonal line s, consisting of vertices and edges of $, 
with the property that it meets every vertical line at exactly one point. Since a extends from 
x = ~°° to x = +oo, any element of the subdivision that is not part of a is either above or 
below it. See figure 10. The elements of a have pairwise disjoint projections on the x-axis, 
and so can be ordered from left to right; the first and last elements are infinite edges. 




Figure 10. A separator. 

A complete family of separators for a monotone subdivision S with n regions is a sequence 
of n - 1 distinct separators Si « s 2 « • • • « s n _ x . There must be at least one region 
between any two consecutive separators, and also one below ai and one above s„-i- If a 
region is below a separator s i} it must also be below any separator « 3 with j > i. We can 
conclude that, if S admits a complete family of separators, its regions can be enumerated as 
Ro, Ri,..., R n -i in such a way that R, « Sj if and only if t < j; in particular, 

Ro « si « i?i « a 2 « ••■ « « n _i « i?„_i. (1) 

Given a complete family of separators and an enumeration of the regions as in (1), let us 
denote by index(fl) the index of a region R in the enumeration. Then s index(K) « R « 
s index(«)+i (whenever those separators are defined). Also, if we let above(c) be the region 
above the edge or vertex e, and below(e) the one below it, the following holds: 

Lemma 7. If i = index(below(e)) and j = index(above(e)) ; then the separators containing 
e will be exactly s i+l ,s i+2 , • ■ ■ , 3j. I 
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6. Existence of a complete family of separators 

It is a well-known result that any acyclic relation over a finite set can be extended to a linear 
(that is, total) order. Therefore, by using lemma 5 we conclude that there is an enumeration 
Ro, . . . , R n -i of the regions of S that is compatible with «, i.e., such that Ri « Rj im- 
plies i < j. Furthermore, any enumeration Ro , Ri , . . . , i2 n _ i of the regions that is compatible 
with -< is also compatible with «, and vice- versa. 

Since a region with nonempty bottom frontier is immediately above some other region, and 
therefore not minimal under -<, the first region in any such enumeration has no bottom 
frontier, and extends to infinity in the — y direction. Since vertical edges are not allowed, its 
i-projection is the whole x-axis. Similarly, the last region Rn-i is the only one with no top 
frontier, and also projects onto the whole z-axis. Therefore, we always have Ro « Ri « 
Rn-i for 0 < i < n - 1. 

We are ready to prove the main result of this section: 

Theorem 8. Every monotone subdivision S admits a complete family of separators. 

Proof: Let Ro, Ri, . . . , R n -i be a linear ordering of the regions of S that is compatible 
with the « relation, i.e. Ri « Rj only if i < j. For t = 1,2, . . . , n — 1, let be the 
collection of all edges and vertices that are on the frontier between regions with indices 
< t and regions with indices > t. For example, figure 10 shows the separator sg. 

Now consider any vertical line /, and let Ri lt Ri 3 , . . . , R* q be the regions it meets, from 
bottom to top. Since / meets Ro and R n -i, and Ri t « Ri, « . . . « Ri t , we will have 
0 = i i < t'2 < ••• < i q = n — 1. Therefore, there is exactly one point on / that is on 
the frontier between a region with index < t and a region with index > t, that is, on 
Furthermore, the intersection of / with s< will be equal to or above that with s,_i, and for 
some such lines (those that meet Ri-i) the two intersections will be distinct. So, we have 
«i « «j « • • • « s„_i. 

Clearly, the elements of Si have disjoint z-projections, and therefore can be ordered from 
left to right; they must be alternately edges and vertices, the first and last being infinite 
edges. To prove that s^ is a separator, it remains only to show that Si is connected; if that 
were not the case, we would have some vertex v of s< that is distinct from, but has the 
same z-coordinate as, one endpoint of an adjacent edge e of s^; see figure 11. 




Figure 11. 

But then we would have u « R « v (or vice-versa), for some region R; and therefore 
e « R « v (or vice- versa), contradicting the construction of s t . Therefore, each s< is a 
separator, and «x > «2 > • • • , *n-i is a complete family of them. I 
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7. A point location algorithm 

Later on we will tackle the problem of efficiently computing a complete family of separators 
for a monotone subdivision. Let's therefore assume for now that we have such a family 
»i, «2, • • • , with a corresponding enumeration of the regions Rq, Ri,..., R n -i satisfying 
(1); we will show next how they can be used to determine, in time 0(log 2 m), the unique 
element of S that contains a given point p. 

The algorithm we will describe is essentially that of Lee and Preparata [LP], and uses two 
levels of binary search. The inner loop takes a separator 8i (as a linear array of edges and 
vertices, sorted by z-coordinate), and determines by binary search an edge or vertex e of s, 
whose x-projection contains the abscissa p x of p. By testing p against e, we will know whether 
p is above a, or below s< (or, possibly, on e itself, in which case the search terminates). The 
outer loop performs binary search on i, so as to locate p between two consecutive separators 
Si and 8i + i, that is to say in a region Ri. 

Besides the separators, the algorithm assumes the index index(i?) can be obtained for any 
given region R, and similarly the adjacent regions above(e) and below(e) can be obtained 
from an edge e, all in constant time. We will see that the construction of these tables is part 
of the process of constructing the family of separators. The search algorithm uses these tables 
to reduce substantially the storage requirements (and also speed up the search a little bit). 
Let T be the infinite, complete binary search tree with internal nodes 1,2,3, .. . and leaves 
0, 1, 2,3,..., as in figure 12. The tree T is used as a flowchart for the outer loop of the search 
algorithm, with the convention that each internal node » represents a test of p against the 
separator s, , and each leaf j represents the output "p is in Rf . While reading the algorithm 
it should be borne in mind that "left" in the tree corresponds to "down" in the subdivision. 
The left and right children of an internal node k of T will be denoted by l(k) and r(fc), 
respectively. We let lca(t',j) be the lowest common ancestor in T of the leaves t and j, that 
is, the root of the smallest subtree of T that contains both » and j. 




Figure 12. The tree T. 



When testing p against a separator, we adopt the convention that each edge contains its 
right endpoint but not its left. This is unambiguous since there are no vertical edges. If the 



Point Location 



12 



algorithm detects that p lies on some edge e during a discrimination against a separator, it 
can terminate the search and, by comparing p with the right endpoint of e, determine if our 
point is a vertex of the subdivision. 

Algorithm 1. Point location in a monotone subdivision. 

{ The algorithm returns in the variable he a reference to the vertex, edge, or region of S containing p.} 

1. Set t +- 0, j <- n - 1, Jfc <- lca(0,n - 1). 
2. While i < j, do: 

{At thit point p is above the separator Si and below ay+i (whenever those separators exist). That is, 
p is in one of the regions Rf, iij+i , . . . , Rj (or in some edge or vertex between two of these regions). 
At each iteration of this loop, either i increases or j decreases, and k (a common ancestor of i and j) 
moves down one level in the tree T.} 

3. If i < k < j then 

4. Find (by binary search) an edge e in Sk such that p x £ lie. 
Let a <— index(above(e)), b «— index(below(e)). 

{By lemma 7, e belongs to the separators 'b+if'k+a «<»• Therefore, by testing p against e 

we conclude it is either on e, above s a , or below «6+l-} 

5. If p is on e, set loc «— e and terminate the search. 

6. If p is above e, set t *— a; else set j *— b. 
7. Else 

[8. If A: > j set k <- /(jfc); else (if jfc < t) set k *- r(Jfc). 
9. Set loc <— Ri and terminate the search. 

The binary search along each separator a* can be performed in O(logm) time if the edges 
of «k are stored in a linear array or balanced binary search tree sorted from left to right. 
By the first iteration, the variable k = lca(0, n — 1) points to an internal node of T at level 
[logn]; at each iteration it descends one level, so we have O(logn) iterations, and a total 
time bound of O (logn log m) = 0(log 2 m). This bound is tight: figure 13 shows a subdivision 
with m edges that realizes it. This example has \/m+ 1 regions and a family of y/m disjoint 
separators with y/m edges each. 

Note that by keeping track of the variables t and j we are sometimes able to skip the binary 
search for p x in some separators. This optimization may improve the average running time of 
the algorithm in practice, but does not affect the worst-case bound. It was included primarily 
to make algorithm 1 more similar to the variants developed further on. 

If we were to independently represent each separator as a linear array, with all its edges and 
vertices, we would have to store n — 1 separators, whose average length can be as large as 
0(m) for some classes of subdivisions. So, the storage requirement of this basic method is 
0(m 2 ) in the worst case. However, after p has been tested against the edge e in step 6, t and 
j are updated in such a way that we will never again look at the edges of any other separator 
that contains e. Therefore, an edge need only be stored in the first separator containing it 
that would be encountered in a search down the tree T. Specifically, if the edge e is in the 
common frontier of regions Ri (below) and Rj (above) , by lemma 7 it belongs to separators 
...Sj and so it suffices to store it in a*, where k is the least common ancestor of t and 
j. This is the highest node in T whose separator contains e. 

Note that only those edges assigned to according to the above rule are actually stored in 
such a structure. In general these will form a proper subset of all the original edges of Sfc, 
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Figure 13. A subdivision that is bad for algorithm 1. 

so between successive stored edges of a fc there may be gaps. See figure 14. Actually, it may 
happen that all the edges of a* are stored higher up in the tree, so that 8k is reduced to a 
single gap, extending from x = -oo to x = +oo. The ordered list of stored edges and gaps 
corresponding to separator will be termed the chain c*. Note that to each subtree of T 
rooted at a node k there corresponds a "partial subdivision" consisting of a range of regions 
and the edges and vertices between them. The separator a* splits this partial subdivision 
into two others, each having half the regions; the gaps of Ck correspond to edges of «*, that 
lie on the upper or lower boundary of the partial subdivision, as shown in figure 14. 




Figure 14. Stored edges and gaps. 

The total storage required to represent the chains is only O(m); in section 10 we show how 
they can be constructed in O(m) time. The derivation of this bound is contingent on our 
ability to compute the least common ancestor of any two leaves of T in 0(1) time. This 
is made possible by the fixed, regular structure of the search tree T (see section 11). The 
point location phase proper could easily be adapted to search a more conventional linked tree 
structure, but such structures seem to admit no simple O(l) algorithm for lea determination. 



Point Location 



14 



8. A faster point location method 

Our hope to obtain a faster algorithm for point location comes from the fact that there is 
some obvious information loss in the method of section 7. Specifically, when we discriminate 
a point against a chain c k , we must localize it in the x coordinate to within an edge or a 
gap of c fc . Yet when we continue the search down some child of k, we start this localization 
process all over again. It would be nice if each edge or gap in a chain pointed to the place on 
each child chain where the x search on that child will finish. The trouble is that an edge or 
gap of the parent can cover many edges or gaps of the child. 

The novel idea in the technique we present in this section is to refine the chains so that 
the localization of p x in one chain allows us to do the same localization in its children with 
only constant extra effort. In its ultimate form, this idea leads to breaking up each chain 
at the x coordinates of all vertices of the subdivision. A bit of thought will convince the 
reader, however, that such an extensive refinement will require quadratic storage for its 
representation. Instead, we describe below a refinement that produces for each chain c k a 
list L k of x-values, defining a partitioning of the x axis into x-intervaU. Each such interval 
of L k overlaps the x-projection of exactly one edge or gap of c fc and at moat two x-intervals 
of the lists Z-t(fc) and £,•(*)• As we will see in section 9, this last condition is compatible with 
keeping the overall storage linear. 

The lists L k and their interconnections can be conveniently represented by a linked data 
structure that we call the layered dag. This is a directed acyclic graph whose nodes correspond 
to tests of three kinds: x-tests, edge testa, and gap teats. A list L k is represented in the dag 
by a collection of such nodes: each x-value of L k gives rise to an x-test, and each interval 
between successive x-values to an edge or gap test. See figure 15. 



gap edge gap edge edge gap 
— > < > < > < >< — > < 




Figure 15. The dag nodes for list L k . 
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An x-test node t contains the corresponding x- value of L k} denoted by xval(t), and two 
pointers left(t) and right(t) to the adjacent edge or gap nodes of L k . An edge or gap test node 
t contains two links down(i) and up(t) to appropriate nodes of L I(fc) and L r(fc) . In addition, 
an edge test contains a reference edge(t) to the edge of c k whose projection covers the x- 
interval represented by t. A gap test node contains instead the chain number chain(t) = k. 
The various types of nodes present are illustrated in figure 16. 




Figure 16. The nodes of the layered dag. 

Let us define more precisely the meaning of the links down(t) and up(*). The properties of 
the refined lists ensure that the x-interval of L k corresponding to an edge or gap test t covers 
either one x-interval / of Lj(fc), or two such intervals I u 7 2 separated by some element x k of 
L|( fc ). In the first case, down(*) points to the edge or gap test of L|( fc) corresponding to the 
interval /; in the second case, down(t) points to the x-test corresponding to the separating 
abscissa x k . Similarly, the link up(t) points to a node of L r(fc ) defined in an analogous manner. 
In the special case when r(Jfc) and l{k) are leaves of T, we let down(f) = up(t) = nil. 
The layered dag then consists of the test nodes for all lists L u L 2 , . . . , £„_!, linked together 
in this fashion. We can use this dag to simulate algorithm 1; the main difference is that each 
time we move from a separator to one of its children, the down or up pointers (possibly 
together with an x-test) allow us to locate x p in the new chain in O(l) time. As before, the 
variables i and j keep track of the interval of separators in which the point p is known to lie, 
namely above Si and below s J+1 . This interval is updated after each edge test exactly as in 
the previous algorithm, and is used during a gap test. When we come to a gap test we know 
that the chain number of the gap is not interior to the interval in question. This gives us an 
unambiguous test as to whether we are above or below the chain of the gap. By the time the 
search algorithm gets to a null up or down link, the interval will have been narrowed down 
to a single region. 

We have now presented enough details about the structure of the layered dag that we can 
give the code for the point-location algorithm. The layered dag contains a distinguished node 
root where the point location search begins. This node is the root of a balanced tree of x-tests 
whose leaves are the edge tests corresponding to the list for the root node of T. 
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Algorithm 2. Fast point location in a monotone subdivision. 

{ This algorithm takes as input the root node of the layered dag and the number n of regions. Its output, as 
in algorithm 1, is placed in the variable loe.} 

1. Set t <— 0, j <— n - 1, t *- root. 

2. While t < j do: 

{At this point we know p is above the separator a,- and below (whenever those separators exist). 

That is, p is in one of the regions Ji,-, JZ,- +1 , . . . , .R^ (or in some edge or vertex between two of these 
regions). The variable t points to a test node in the layered dag, which together with its descendants will 
allow us to locate the point p among those regions.} 

3. If £ is an edge test then let e <— edge(t) and do: 

{At this point we know p x lies within lie.} 

4. If p is on e, set loc <— e and terminate the algorithm. 

5. If p is above e, set t <— up(t) and t <— index(above(e)). 
Else set t +— down(t) and j «— index(below(e)). 

6. Else if t is an x-test then do: 

{ The following x-test routes us to the appropriate edge of the next chain we need to test against.} 
7. If p x < xval(f) then t <- left(t) else t <- right(i). 
8. Else t is a gap test; do 

{ We have already compared p against the appropriate edge of the chain of the gap test We just 
need to reconstruct how that comparison went} 

9. If j < chain(t) then t *— down(f) else t <— up(t). 
10. Set loc «— Ri and terminate the search. 



9. Computing the layered dag 

Now that we understand how the layered dag is to be used, we will describe how it is to be 
constructed. Our starting point will be the tree T and the chains Ck defined in section 7; 
recall that Ck consists of those edges of Sk that do not belong to any ancestor of Sk (that is, 
to any separator whose index is an ancestor of k in T). 

Our construction of the layered dag proceeds from the bottom up and happens simultaneously 
with the refinement of the chains c^- We first describe how the x values in Lk are obtained. 
Note that we already have at our disposal three sorted lists of i values: those corresponding 
to £j(fc), to ^r(fc)) an d also to the endpoints of the edges stored with the chain cjt associated 
with node k in the chain tree. The x values in Lk are a merge of those in cjt, and every other 
one of those present in each of L^k) and -£<r(Jfc)- By convention, if A: is a terminal node of the 
chain tree (so it corresponds to a region), or A: > n, then Lk is empty. We imagine now that, 
in a bottom-up fashion, this is done for every node k of the chain tree. The propagation of 
every other value from the children to the father constitutes the chain refinement we had 
mentioned in section 8. This refinement has two desirable properties: 
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Lemma 9. An interval between two successive x values in Lk overlaps at most two intervals 
,n L i{k), or i r (fc), respectively. | 

Lemma 10. The total number of x values in the lists Lk, summed over all k, is at most 
4m. 

Proof: If ak denotes the number of edges in Ck, then 

X) a * = m > 

fcGT 

since each edge of the subdivision occurs in exactly one chain c*. Let 6* denote the number 
of x values in L k , and A k (resp. B k ) denote the sum of a< (resp. &,) over all nodes i in the 
subtree rooted at A;. To prove the lemma it suffices to show that 

B r < 4A r = 4m 

for the root node r = lca(0,n - 1). 

As an inductive hypothesis now assume that 

Bi + bi < 4Ai (2) 

for i = l(k) or t = r(k). This hypothesis is trivially true for the leaves of T. Observe now 
that 

Bk = 5|(fc) + B r (fc) + b k , 

and that 

b k < 2a fc + {b i(k) + K {k) )/2, 

since each edge in c k contributes at most two x- values to Lk- Applying the inductive 
hypothesis (2) yields 

B k + b k < 4A fc) 

which proves the same conclusion for k. This concludes the proof of the lemma. | 

Intuitively, half of the x values of c k propagate to the father chain, a quarter to the grandfa- 
ther, an eighth to the grand-grandfather, and so on. Although this argument is not strictly 
correct, it illustrates why we can expect the combined length of the L lists to remain linear 
in m, and in fact just twice as much as the total number of x values in the chains Cfc. 

The construction of the x-test, edge test, and gap test nodes representing the list Lk in the 
layered dag is now straightforward, as well as the setting of the left and right fields of the x- 
tests. The down links of L fc can be set according to the rules stated in section 8, by traversing 
simultaneously the two lists Lk and from left to right. The up pointers are analogous. 
In fact, it is possible to build Lk and link it to the rest of the dag in a single simultaneous 
traversal of Cfc, £j(fc), and £ r (fc)- This bottom-up process terminates with the construction of 
the root chain L r , where r = lca(0, n - 1). As a final step we produce a tree of i-test nodes 
corresponding in a natural way to a binary search on the a:- values of the list L r . The leaves of 
the tree are made to point to the appropriate edge test nodes of L r . All nodes of the layered 
dag can be reached from the root of that tree. 
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Algorithm 3. Construction of the layered dag. 
1. Set t «- 1. While t < n do: 

{Construct one more level of the tree T. The first node in this level is i and the difference between 
successive nodes is 2t. A node k on this level is the common ancestor of the leaves k - i through 
k + i- 1; its children (ifi> \) are the internal nodes k - i/2 and k + i/2.} 

2. Set A; <- t\ While it - i < n do: 

{Create the list L k from the chain c fc and the lists V = L, (fc) and L" = Z, r(fc) (if they exist).} 

3. If t = 1, set V *— <f>, else set V <— L k _i/ 2 - 

4. If t = 1 or k + i/2 > n, set L" <- <f>, else set L" *- L k+i/2 . 

5. If k > n, let c k be a single gap from i -= — oo to x = +oo. 

Split the edges and gaps of c k at every other x- value of V and L". 

6. Set L k <— <f>. For each edge or gap e in c k , do: 

7. Append to L k an edge or gap test t representing e. 
Set edge(f) «— e or chain(t) «— A:, as appropriate. 

8. If V = <t>, let down(f) «- n»7. 

Else, if lie overlaps only one x-interval of L', let down(t) point to the 
corresponding edge or gap test of V . 

Else lie overlaps two i-intervals of L'; create a new z-test node t' that 
chooses between the two corresponding edge or gap tests of L', and let 
down(f) «- t'. 

9. Similarly, set up(t) to nil, to an edge or gap test of L", or to a new 
x-test that chooses between two tests of L" . 

10. Set k <- k + 2i. 
11. Set t <- 2i. 

12. Let r «— lca(0,n - 1). Construct a binary tree of x-tests for the list L r , and 
let root point to this tree. 

Note that several different nodes of L k may point to the same node of a child; an example 
is shown in figure 17. Thus the resulting dags are not trees. We remark that the structure 
built by Kirkpatrick [Ki] also corresponds to a dag. This "sharing" of subtrees seems to be 
an essential feature of algorithms for point location that simultaneously attain optimal space 
and query time. 




Figure 17. convergence in the dag. 



Point Location 



19 



To cut down on the number of links, we may consider storing the edge and gap test nodes 
of each list Lk in consecutive memory locations, in their natural left-to-right order, with 
the i-tests between them. The initial location of x p in L r could then be carried out by 
standard binary search. Also, in the construction of Lk we would be able to scan the lists 
Z/j( fc ) and L r ( fc ) without any auxiliary pointers or tables. Finally, this sequential allocation 
would eliminate the need for the left and right links of x-tests, a fact that may reduce the 
total storage used by pointers in the dag by about one third. 

In any case, the initial location of x p in the root list L r can be determined in 0(log m) time. 
After that, algorithm 2 executes exactly [lg n] edge or gap tests (one at each level of T), and 
at most that many x-tests. So the total query time is O(logm). A list Lk with t x- values is 
represented in the layered dag by at most t x-tests and t + 1 edge/gap tests and, as we have 
seen, it can be constructed in 0{t) time. Using lemma 10 we conclude that the layered dag 
contains at most 4m x-tests and 4m + n — 1 edge and gap tests, and can be built in total 
time 0(m + n) from the chain tree T. In summary, we have shown that 

Theorem 11. Assuming that the chain tree for a subdivision with m edges is given, the 
layered dag data structure can be constructed in 0{m) time, takes O(m) space, and allows 
point location to be done in O(logm) time. I 
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10. Constructing a complete family of separators 

We proceed now to the description of an efficient algorithm that constructs the chain tree T 
representing a complete family of separators for a given monotone subdivision S . As suggested 
by the proof of theorem 8, the first pass of this algorithm enumerates the regions in a linear 
order compatible with «. A second pass enumerates the edges and vertices of S, in such 
a way that the edges and vertices of each separator are visited in order of increasing in- 
coordinate (even though this may not be true for elements belonging to different separators). 
Therefore, by just appending each visited element to its appropriate separator, we will get 
the sorted lists required in algorithms 1 and 3. 

As in section 3, we will add to S two dummy vertices at x = — oo and x = +oo, which are 
endpoints of all edges with infinite left or right extent, respectively. If we orient every edge 
of S from right to left, we obtain a planar embedding of a directed, acyclic graph S with 
exactly one source and one sink (solid lines in figure 18). The enumeration we need in the 
second pass is basically a compatible traversal of this graph, in which we visit a vertex only 
after visiting all its outgoing edges, and we visit an edge only after visiting its destination. 
This is a form of topological sorting, as discussed in [Kn] . 




Figure 18. The graphs 5 and 5*. 

Consider now the dual graph S* whose vertices are the regions of S, and where for each edge 
e of S there is an edge c' from above(e) to below(e). By what we have seen in sections 4-6, S* 
too is acyclic and has exactly one source and one sink. We can draw S* on top of S (dotted 
lines in figure 18) in such a way that each of its edges e' crosses only the corresponding edge 
e of S, and that exactly once (going down). Therefore, 5* is planar, and corresponds to 
the topological dual of S. It turns out that S* represents for the first pass what S does for 
the second: a compatible traversal of S* will visit the regions in an order consistent with -< 
and «. 

Therefore, both passes reduce to a generic graph traversal algorithm, applied first to S* and 
then to S. In the first pass, as each region R is visited, we store in a table (or in a field of the 
region's record) its serial number index(i2) in the enumeration. In the second pass, as each 
edge or vertex e is visited, we obtain the indices a = index(above(e)) and b = index(below(e)) 
of the regions immediately above and below it, and append e to the separating chain c k where 
k = lca(a, b). With an appropriate representation for subdivisions, it is possible to accomplish 
this graph traversal without any auxiliary storage. This topic is discussed in section 12. 
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11. The lowest common ancestor function 

If the construction of the separating chains c, as described above is to run in 0(m) time, it 
is essential that the total time to compute lca(index(below(e)),index(above(e))) for all edges 
e be 0(m). This rules out the straightforward algorithm that starts at the leaves t and j, 
and moves up one level of T at a time, following "parent" links, until the two paths join at 
a common node. This naive algorithm has running time f2(log |t - and it is possible to 
have subdivisions in which |index(below(e)) - index(above(e))| = n(v^n) for f2(m) edges e, 
thus giving an overall running time of fi(m log m). 

Algorithms for computing in O(l) time the least common ancestor function on general binary 
trees have been published by Harel [H]. His algorithms are probably too complex to be of 
practical use here, but they can be considerably simplified thanks to the regular structure 
of our search tree T. On this tree, the value of lca(t,j) has a simple interpretation in terms 
of the binary representations of t and j. Let v = [lg n] be the number of bits needed to 
represent any number from 0 through n - 1, and let 

x = (a 1/ _ 1 a„_ 2 ... o r+2 0 o! r ... a\ a' 0 ) 2 , 

j = (a„_! a„_ 2 . . . a r+2 1 a' r ' . . . a'{ a#) 2 . 

Then 

lca(i,j) = (a„_! a r+2 1 0 ... 0 0 ) 2 . 

Loosely speaking, lca(i',jf) is the longest common prefix of j and j, followed by 1 and padded 
with zeros. 

An efficient formula for computing lca(t',j) is based on the function msb(fc) = lca(0,fc) = 
2 pg *1 -l ) tQ g mQ8t s ig n ifi cant ^ 0 j ^ i ta va i ues for A; = 1, 2, . . . , n - 1 are 1, 2, 2, 4, 4, 4, 4, 
8, 8, . . . ,2" -1 ; these numbers can be easily precomputed in O(n) time and stored in a table 
with n - 1 entries, so we can compute msb(fc) in O(l) time. Then we can express the lea 
function as 

lca(i, j) = j A -.(msb(t © j) - 1) 

where ©, A, and -. are the boolean operations of bitwise "exclusive or", "and", and comple- 
ment. We assume these boolean operations can be computed in O(l) time, like addition and 
subtraction. We feel their inclusion in the model is justified, since their theoretical complexity 
is no greater than that of addition, and most computers have such instructions.t For use in 
this formula, it is preferable to tabulate -.(msb(Jfc) - 1) instead of msb(fc). 
Another way of computing lea is based on the bit-reversal function rev(ifc), that reverses 
the order of the last v bits in the binary expansion of k. For example, for r» = 16 and 
A: = 0, 1, . . . , 15, the values of rev(fc) are 0, 8,4, 12, 2, 10,6, 14, 1, 9, 5, 13, 3, 11, 7, 15. Using this 
function we get the formula 

lca(t, j) = rev(fc © (k - 1)) A j, 

where k = rev(t © j). 



f ^ t >. 8ome . machlnes ca n compute \\gk] - 1 from k ("find first 1 bit") and 2« from q ("shift left q 
bits ), m a single instruction cycle. Under a more rigorous log-cost complexity model, all time and space 
bounds m this paper and in the main references should be multiplied by an extra logm factor. 



Point Location 



22 



12. Compatible traversal of an acyclic graph 

The input to the compatible graph traversal routine we mentioned in section 10 is a planar 
embedding G of a directed, acyclic, and connected graph with exactly one sink and one 
source, both on the exterior face of G. See figure 19. 




sink 



Figure 19. The input to algorithm 4. 

Such an embedding defines a "counterclockwise" circular ordering of the edges incident to 
any given vertex u. The post-order traversal of G is defined as a listing of all its vertices and 
edges such that 

(i) an edge is listed (or visited) only after its destination, 

(ii) a vertex is visited only after all its outgoing edges, and 

(iii) edges with same origin are visited in counterclockwise order. 

This is clearly a compatible traversal as defined in section 10. The post-order traversal is 
unique, and is a particular case of the general depth- first graph traversal described by Tarjan 
[Ta]. This problem admits a straightforward recursive solution. Given a vertex u (initially 
the source of G), we enumerate the edges out of u, in counterclockwise order. For each edge 
e, we first recursively apply the procedure to its destination v (unless it has been previously 
visited). We then visit c and proceed to the next edge. After all edges out of u have been 
visited, we visit u and exit. This algorithm runs in 0(m) time and requires only O(m) 
auxiliary storage, the latter consisting of the recursion stack and one mark bit per vertex (to 
distinguish the nodes that have already been visited). 

In the rest of this section we will show that this post-order traversal can be performed 
without an auxiliary stack or any mark bits on the vertices, provided the data structure 
used to represent the subdivision is rich enough. This improvement is of significant practical 
interest, even though it does not affect the 0(m) space bound. 

As we observed, the embedding of G in the plane defines a counterclockwise ordering of the 
edges incident to a given vertex u. In this ordering all outgoing edges occur in consecutive 
positions, and the same is true of the incoming ones. To see why, consider any two edges 
ei,ei entering u, and any two edges gi,g2 leaving u. Let jti, jt 2 be two paths from the source 
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Figure 20. 

vertex of G that end with the edges ci and e 2 , and let ai,<7 2 be the two paths to the sink 
vertex that begin with <7i,fl 2 . See figure 20. 

Since G is acyclic, both t x and jt 2 are disjoint from a x and a 2 (except for u itself). Now, the 
paths TTx and ?r 2 together divide the plane in two (or more) regions. If the two pairs of edges 
were interleaved, at least one of the paths o x and <r 2 would have to cross K\ or * 2 , since they 
start on different regions but have a common destination. This proves the above assertion. 

If u has both incoming and outgoing edges, this result establishes a linear order for each class 
with well-defined "first" and "last" elements, which will be denoted by first in(u), lastin(u), 
first out(u), and lastout(u). See figure 21. To make this definition meaningful also for the 
source and sink of G, we will introduce a dummy edge base(G) that connects the sink back to 
the source across the exterior face of G. We may consider the resulting graph G' as embedded 
on the upper half of a sphere, with base(G) running across the bottom half, as in figure 22. 
In the graph G', the first and last outgoing edges of the source of G will be those incident to 
the exterior face of G. A similar statement applies to the sink of G. 




Figure 21. 



Figure 22. 
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Let us mechanically translate the recursive post-order algorithm into an iterative one, using 
an explicit stack Q to save the value of e when simulating recursive calls (the value of u need 
not be saved, since it is always the origin of e). 

Algorithm 4. Post-order traversal of an acyclic planar graph with one source and one sink. 



5. If e 7^ base(G), visit e. 

6. While e = last out(u) do 

7. Visit u. 

8. If Q is empty, the algorithm terminates. Else, 

9. Pop the topmost edge of Q, and assign it to e. 



10. Visit e. 
11. Set e to the next edge out of u. 

The following claims about algorithm 4 are a direct consequence of the planarity and acyclic- 
ity of G, and can easily be proved by induction on the number of executed steps: 

(I) Before every step, the edges in the stack Q form a path in G from the source to u; 

(II) an edge is stacked onto Q (step 4) only if its destination is still unvisited; 

(III) an edge is unstacked (step 9) only after its destination has been visited; 

(IV) an edge is visited only after its destination has been visited; 

(V) a vertex is visited only after its last outgoing edge has been visited; 
(VI) every edge is stacked at most once; 

(VII) for any given vertex, at most one incoming edge ever gets stacked; and 

(VIII) an edge is visited only after the previous edge with same origin (if any) is visited. 

In particular, from (IV), (V), and (VIII) we conclude that all vertices and edges of G are 
visited, and conditions (i)-(iii) are satisfied. Algorithm 4 defines a subgraph H of G, consisting 
of all the edges that ever get stacked onto Q. Every vertex of G is reachable from the source 
via a directed path in H (the contents of the Q at any instant when the variable u is that 
vertex), and has at most one incoming edge in H. It follows that H is an oriented spanning 
tree of G, as illustrated in figure 23. We remark that the order in which the vertices of G 
are visited corresponds to the post-order traversal of the tree H, as defined by Knuth [Kn]. 
We will show now that the only edge of H (if any) entering a vertex u is lastin(u). More 
precisely, the following lemma holds: 




Set u to the origin of e. 
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source 



sink 



Figure 23. The spanning tree H. 

Lemma 12. Before any step of algorithm 4, if the stack Q is not empty, its topmost edge 
is last in(u); otherwise u is the source of G. 

Proof: The second part of the lemma is obvious, since the contents of Q is always a path 
in G from the source to u. Let us then assume Q is not empty. Let ui be the current value 
of u, 7T be the path in Q, and e y be the last edge of jt (i.e., the top of Q). See figure 24. 




sink 



Figure 24. 

Suppose e 2 is an edge distinct from ei but having the same destination ui. From assertions 
(II) and (III) above we conclude that e x and e 2 have yet to be visited. By the time e 2 is 
visited by algorithm 4, the variable u will contain the origin of e 2 ; the contents of Q at that 
time, plus the edge e 2 , will define another directed path n' from the source to tii. Since G 
is acyclic, the path %' neither contains nor is contained in ic; in fact, because of property 
(VII) the paths tt and tt' must diverge exactly once, at some vertex w ^ ui, and converge 
exactly once at u x . Therefore all edges of jt between w and u x must be unstacked (and 
visited) before e 2 is visited. In particular, the edge a through which jt leaves w is visited 
before the corresponding edge a' of jt', and therefore a' must follow a in counterclockwise 
order around w. 

Let R be the complement of the outer face of G, and let a be any directed path from ui to 
the sink of G. The concatenation of jt and a divides R in two (not necessarily connected) 
regions, which we call left and right (with the direction of xo~ being taken as "forwards"). 
The path jt' cannot cross a (otherwise the two would give rise to a cycle), and its first 
edge a' lies in the left region; therefore, after leaving w the path *' must lie entirely 
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in the left region, and in particular ei precedes e\ in the counterclockwise order around 
ui. By repeating this argument for all possible edges t<i ^ e\ into ui, we conclude that 
ei = lastin(u). | 



Therefore, instead of retrieving e from the stack in step 9, we can simply set e +— last in(tz). 
Note also that the test "e ^ base(G) and v has not been visited yet" in step 3 is evaluated 
and yields true if and only if e is pushed into Q by the following step, so we can replace 
that test by the condition "e ^ base(G) and e = last in(v)" . These observations enable us to 
dispense with the recursion stack and the "visited" bits on the vertices. 

A representation for the embedded graph G that allows the efficient computation of first out 
and its companions is the quad-edge data structure [GS]. This representation is similar to 
the well-known "winged edge" and DCEL data structures [B, MP], but has over them the 
important advantage of encoding simultaneously both G and its dual embedding, in precisely 
the same format, at a negligible extra cost in storage. This allows the post-order traversals of 
both S and 5* to be performed by the same procedure, applied to the same data structure. 

The quad-edge structure by itself can only represent an undirected embedded graph, such 
as the undirected subdivision S. However, every edge e of 5 is represented by two distinct 
records in the structure, corresponding to the two possible orientations of e. Therefore, to 
refer to the edge e of the structure we must actually refer to a specific directed version of e. 
Given such a directed edge c, the quad-edge data structure gives immediate access to: 
org(e) the origin vertex of e 
dest(e) the destination vertex of e 
onext(e) the next counterclockwise directed edge with the same origin, 
dnext(e) the next counterclockwise directed edge with the same destination, 
sym(e) the same edge directed the other way, and 
rot(e) the dual of the edge e, directed from the right to the left faces of e. 

A directed graph G, such as S or S*, can be represented by the quad-edge encoding of the 
corresponding undirected graph, plus a predicate forward(e, G) that tells whether the directed 
edge e of the structure is actually an edge of G. Clearly, e is in G if and only if sym(e) is not 
in G, so forward(e,G) = -> forward(sym(e),G). In our case, forward(e, S) is simply the test 
of whether the x-coordinate of dest(e) is smaller than that of org(e). Similarly, in the dual 
graph S* the predicate forward(e, S*) tests whether the region dest(e) is immediately below 
the region org(e). This turns out to be the same as -iforward(rot(e),S). The dummy edges 
that we must add to S and S* are the only exception: we have rot(base(5*)) = base(S), and 
yet both are forward. For both graphs, we also have 

e = last in(dest(e)) •<=> forward(e, G) A -> forward(dnext(e), G), 
e = last out (org (c)) o forward(e, G) A forward(onext(e), G), 
c = last in(u) => sym(dnext(e)) = first out(u). 



These identities allow us to remove also the calls to first out and its companions from algo- 
rithm 4. Algorithm 5 below incorporates all these modifications. 
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Algorithm 5. Post-order traversal of an acyclic planar graph with one source and one sink 
using O(l) auxiliary storage. 

l.Set u <- dest(base(G)) and e <- sym(dnext(base(G))). 
2. Repeat the following steps: 

{At this point u it unvisited, forward(e, G) it true, and e it the first unvitited edge out of u.} 

3. While e / base(G) A -> forward (dnext(e),G), do 

{Here e it the last edge into its destination v, to v hat not been visited yet} 
4. Set u <— dest(e), and set e <— sym(dnext(e)). 
{Thit sets e to first out(u).} 
{Now e it the dummy edge, or its destination hat already been visited.} 

5. If e ± base(G), visit e. 

6. While ->forward(onext(e),G) do 

{Here e hat been visited and is last out(u).} 

7. Visit U. 

8. If u = dest(base(G)), the algorithm terminates. Else, 

{ Compute last in(u), and backtrack through it.} 

9. Set e *— sym(onext(e)). 

While forward(dnext(e),G), do e «- dnext(c). 
Set u <— org(e). 

10. Visit e. 

11. Set c *— onext(e). 

The equivalence between algorithms 4 and 5 is straightforward. It is also easy to show that 
the latter runs in 0(m) time for a subdivision with m edges. Every edge of the graph is 
assigned to e at most twice: once in the enumeration of the edges out of a vertex t; (step 11), 
and once in the determination of lastin(v) (step 9). 
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13. Conclusions and applications 

We have introduced a new data structure, the layered dag, which solves the point location 
problem for a monotone subdivision of the plane in optimal time and space. The main idea 
has been to refine the chains introduced by Lee and Preparata and connect the refined chains 
by links. The latter concept originates with Willard [W] and has found frequent application 
since then. The layered dag can be built from standard subdivision representations in linear 
time, as follows. We use the graph traversal algorithm of section 12 to enumerate the regions 
of the subdivision in a way compatible with the vertical ordering presented in section 4. 
Another traversal of the subdivision then allows us to build the chain tree representing a 
complete family of separators, as in section 10. Finally, the layered dag is built from the 
chain tree, as explained in section 9. The point location algorithm using this structure has 
been given in section 8. Compared to previous optimal solutions, the advantage of the layered 
dag is that 

• it admits a simple, practical implementation, and 

• it can be extended to subdivisions with curved edges. 

We will not discuss in detail here how to generalize our method to work for curved-edge 
subdivisions. Certain requirements for such a generalization to work are clear. We must be 
able to break up edges into monotone pieces, to introduce the additional edges required by 
regularization, and to test on what side of (a monotone segment of) an edge a point lies. Our 
time bounds will be maintained as long as we are able to in constant time: 

• cut an edge into monotone pieces, 

• add a monotone regularization edge between two existing monotone edges, and 

• test if a point p is above or below a monotone edge e. 

The layered dag also yields improved solutions for several other problems in computational 
geometry. All these problems are reduced to the subdivision search problem treated earlier. 
For example, subdivisions with circular edges occur in the weighted Voronoi diagram of a 
point set [AE]. There, each point p in a finite set U has associated a positive weight w(p) and 
the region R(j>) = {x\ d{x,p)/w(p) < d(x,q)/w(q), for all q G U }. The layered dag offers the 
first optimal method for locating a point in the diagram defined by these regions. 

Finally, certain problems related to windowing a two-dimensional picture given as a collection 
of line segments have been reduced to subdivision search by Edelsbrunner, Overmars, and 
Seidel [EO]. The layered dag provides a way to extend their methods to more general curves 
without losing efficiency. 
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Section 


Meaning 
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preprocessing time, storage, and query time of an algorithm 


n r 
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projection of R on the z-axis 


Ay> B 
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set A is above set B 
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5 


4 


a monotone subdivision 
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an enumeration of the regions compatible with 


index(iZ) 


5 


mucA ui region n in ine compaiiuie enumeration 


below(e), above(c) 


5 


regions adjacent co ine edge e 




7 


** xvwiufc ku uc lULabcd) a ii ci us x -coordinate 


T 


7 
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7 
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/oc 


7, 8 


the element of S containing p 


=fc 


7, 10 


the edges and traDS of sir that do not occur hicrher nn in *P 




8, 9 


list of edge tests, gap tests, and z values representing the chain 


xval(t) 


8 


the z-value stored in the z-test node t 


left(t),right(t) 


8 


pointers in the z-test t to the adjacent edge or gap tests. 


down(t),up(t) 


8 


test nodes in children chains that overlap the edge or gap test t 


edge(t) 


8 


edge of the subdivision represented by the edge test t 


chain(t) 


8 


index of the chain the gap test t belongs Rto 




8 


> * 
a pointer to no record 


root 


8 


root of the layered dag 


msb(A;) 


11 


most significant bit of k 


rev(fc) 


11 


the result of reversing the binary representation of k 


S 


10 


directed graph obtained by oriented all edges of S to the right 


S* 


10 


dual of 5, with all edges directed downwards 


G 


12 


the graph (either 5 or S*) that is to be traversed in post-order 


first in(u) 


12 


first edge (in ccw order) entering vertex u of G' 


last in(u) 


12 


last edge entering u 


first out(u) 


12 


first edge (in ccw order) leaving u 


last out(u) 


12 


last edge leaving u 


base(G) 


12 


a dummy edge connecting the sink of G to its source 


G' 


12 


the graph G plus the dummy edge base(G) 


Q 


12 


a stack used in algorithm 4 


org(e),dest(e) 


12 


origin and destination of the edge E 


sym(e) 


12 


the edge e taken with opposite orientation 


rot(e) 


12 


dual of edge e, oriented from right to left of « 


onext(e) 


12 


edge following e in ccw order around org(e), oriented away from org(e) 


oprev(e) 


12 


edge preceding e in ccw order around org(e), oriented away from org(e) 


dnext(e) 


12 


edge following e in ccw order around dest(e), oriented towards dest(«) 


forward e 


12 


a predicate that tells whether e or sym(e) is an edge of G' 
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acyclic 8a, 10a, 14b, 20, 22a, 

23a, 24, 25b, 27 
balanced (search) tree 6, 12b, 

15b 

bitwise boolean operations 21b 
chain 2, 13, 14, 15, 16, 17, 18a, 

19, 20, 21a, 28a 
chain tree 2b, 16b, 19, 20a, 28a 
complete family of separators 2b, 

9b, 10, 11a, 20a, 28a 
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2, 28a 
dual graph 20b, 26b 
edge test 14b, 15, 16a. 17b, 18a, 
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Euler's Theorem 3a 
firstout 26 
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generalization 2a, 28a 
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regular 3b, 5a, 13b, 21a 
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sweep(ing) line 5b, 6a 
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Willard 28a 
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